Praca domowa nr 5 dotyczy wyjaśnień globalnych modelu przy użyciu metod odpowiadających na pytanie, jak poszczególne zmienne wpływają na średnią predykcję. Wykorzystane metody to PDP - Partial Dependence Profiles i ALE - Accumulated Local Effects.
PDP pokazują jak zachowuje się wartość oczekiwana predykcji modelu jako funkcja wybranej zmiennej objaśniającej, co można estymować jako średnią z profili CP dla każdej obserwacji. Jednak PDP mogą być mylące w przypadku skorelowanych zmiennych objaśniających - np. rozważając nasz zbiór nierealistyczne jest rozważanie rezerwacji na 0 nocy weekendowych i ponad 5 nocy w dni powszednie.
Dlatego warto wykorzystać ALE, które, wykorzystując rozkład łączny zmiennych objaśniających innych od wybranej, intuicyjnie pokazują jak zachowuje się wartość predykcji w małym oknie wokół wartości wybranej zmiennej. (Ciężko opisać to bez równań :) )
import pandas as pd
import numpy as np
import pickle
import dalex as dx
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import LabelEncoder, OneHotEncoder
from sklearn.compose import ColumnTransformer
from sklearn.impute import SimpleImputer
from sklearn.pipeline import Pipeline
from sklearn.ensemble import GradientBoostingClassifier
from xgboost import XGBClassifier
import warnings
warnings.filterwarnings('ignore')
# Wczytanie danych i wstępne operacje
full_data = pd.read_csv("data/hotel_bookings.csv")
full_data["agent"] = full_data["agent"].astype(str)
treshold = 0.005 * len(full_data)
agents_to_change = full_data['agent'].value_counts()[full_data['agent'].value_counts() < treshold].index
full_data.loc[full_data["agent"].isin(agents_to_change), "agent"] = "other"
countries_to_change = full_data['country'].value_counts()[full_data['country'].value_counts() < treshold].index
full_data.loc[full_data["country"].isin(countries_to_change), "country"] = "other"
# Określenie cech uwzględnionych w modelach
num_features = ["lead_time", "arrival_date_week_number",
"stays_in_weekend_nights", "stays_in_week_nights",
"adults", "previous_cancellations",
"previous_bookings_not_canceled",
"required_car_parking_spaces", "total_of_special_requests",
"adr", "booking_changes"]
cat_features = ["hotel", "market_segment", "country",
"reserved_room_type",
"customer_type", "agent"]
features = num_features + cat_features
# Podział na zmienne wyjaśniające i target
X = full_data.drop(["is_canceled"], axis=1)[features]
y = full_data["is_canceled"]
# Label encoding - transformacja danych
categorical_names = {}
X_le = X.copy()
for feature in cat_features:
col = X_le[[feature]]
cat_transformer = SimpleImputer(strategy="constant", fill_value="Unknown")
col = cat_transformer.fit_transform(col)
X_le[feature] = col
le = LabelEncoder()
le.fit(X_le[[feature]])
X_le[[feature]] = le.transform(X_le[[feature]])
categorical_names[feature] = le.classes_
# Podział zbioru z label encodingiem
X_train_le, X_test_le, y_train, y_test = train_test_split(
X_le, y,
test_size=0.2, random_state=42)
Pod uwagę wezmę główny model, który opracowujemy w ramach grupy - las losowy z label encodingiem. Kod można znaleźć na repozytorium. W późniejszym etapie rozwiązania wczytam i wytrenuję również inne modele dla porównania wyników.
rf_le_model = pickle.load(open("RF_model", "rb"))
rf_le_explainer = dx.Explainer(rf_le_model, X_train_le, y_train, label = "RF LE")
rf_le_pdp = rf_le_explainer.model_profile(N=1000, random_state=42)
rf_le_pdp.plot()
Profile PD dla naszego głównego modelu potwierdzają wcześniejsze przypuszczenia i wnioski odnośnie tego, jak poszczególne zmienne wpływają na predykcję.
Duże zróżnicowanie widać dla zmiennej lead_time. Dla mniejszego czasu poprzedzającego rezerwację prawdopodobieństwo jej odwołania jest mniejsze. Wartość ta rośnie bardzo szybko aż do około miesięcznego okresu lead_time. Dla okresu poprzedzającego o długości powyżej około 300 dni model daje już predykcję odwołania rezerwacji. Zatem można podejrzewać, że nie warto jest umożliwiać rezerwacje z tak dużym wyprzedzeniem.
Widać też, że dłuższe pobyty nieco częściej są odwoływane aniżeli te krótkie. Nie jest to znacząca różnica, jednak zauważalna. (stays_in_week_nights, stays_in_weekend_nights)
Bardzo dobrze widać również, że fakt wcześniejszych odwołań rezerwacji (previous_cancellations) wpływa na predykcję kolejnego odwołania, co również wydaje się naturalne.
W przeciwną stronę działa zmienna oznaczająca liczbę wcześniejszych nieodwołanych rezerwacji (previous_bookings_not_canceled), natomiast niezerowe wartości nie zmieniają wartości predykcji tak bardzo.
Innymi czynnikami, które wpływają na predykcję nieodwołania rezerwacji są podejmowane przez klienta dodatkowe działania w sprawie bookingu - dodatkowe prośby (total_of_special_request), zmiany w samej rezerwacji (booking_changes) i w największym stopniu - rezerwacja miejsca parkingowego (required_car_parking_spaces) (możemy myśleć o tym, że w tych hotelach jest to dodatkowo płatne albo też osoby podróżujące samochodami są mniej zależne od środków komunikacji publicznej, więc ich przyjazd jest pewniejszy).
Profile zostały wygenerowane również dla zmiennych kategorycznych, jednak w przypadku naszego modelu i wykorzystanego label encodingu są one przedstawiane jak wykresy dla zmiennych ciągłych. Natomiast widzimy, że niektóre z nich są mocno zróżnicowane - te, które zidentyfikowaliśmy przy użyciu permutacyjnej ważności zmiennych jako zmienne znaczące. Po raz kolejny widzimy także, jak zmienia się predykcja modelu, jeśli osoba rezerwująca jest z Portugalii.
Profile PDP są przydatne także w porównaniach różnych modeli, dlatego sprawdzimy jak wyglądają one w przypadku lasu losowego z one hot encodingiem (modelu, od którego zaczynaliśmy projekt) oraz modeli XGBoost i GradientBoostingu.
# Podział zbioru do one hot encodingu (wykonywanego w modelu)
X_train_oh, X_test_oh, y_train_oh, y_test_oh = train_test_split(
X, y,
test_size=0.2, random_state=42)
# Definicja transformacji danych w modelu z one hot encodingiem
num_transformer = SimpleImputer(strategy="constant")
cat_transformer = Pipeline(steps=[
("imputer", SimpleImputer(strategy="constant", fill_value="Unknown")),
("onehot", OneHotEncoder(handle_unknown='ignore'))])
preprocessor = ColumnTransformer(transformers=[("num", num_transformer, num_features),
("cat", cat_transformer, cat_features)])
# Import modelu RF
rf_ohe_model = pickle.load(open("rf_model_old_version", "rb"))
xgb_model = XGBClassifier(random_state=42, n_jobs=-1)
xgb_model.fit(X_train_le, y_train)
gb_model = GradientBoostingClassifier(random_state=42)
gb_model.fit(X_train_le, y_train)
rf_ohe_explainer = dx.Explainer(rf_ohe_model, X_train_oh, y_train_oh, label = "RF OHE")
xgb_le_explainer = dx.Explainer(xgb_model, X_train_le, y_train, label = "XGB LE")
gb_le_explainer = dx.Explainer(gb_model, X_train_le, y_train, label = "GB LE")
rf_ohe_pdp = rf_ohe_explainer.model_profile(N=1000, random_state=42)
xgb_le_pdp = xgb_le_explainer.model_profile(N=1000, random_state=42)
gb_le_pdp = gb_le_explainer.model_profile(N=1000, random_state=42)
rf_le_pdp.plot([rf_ohe_pdp, xgb_le_pdp, gb_le_pdp])
Zauważmy, że profile dla lasów losowych z obiema wersjami kodowania zmiennych kategorycznych są niemal identyczne (na wykresach zielona i czerwona linia).
Patrząc szerzej, profile są porównywalne dla wszystkich modeli. Najbardziej znaczące różnice widać w zmiennych dotyczących wcześniejszych rezerwacji danego klienta - ilości odwołanych i nieodwołanych rezerwacji. Modele XGB i GB bardziej "premiują" niezerowe wartości - różnice predykcji są większe. W ogólności to właśnie model XGB jest najbardziej wrażliwy na wszelkie zmienne, co widzimy po najbardziej poszarpanym kształcie krzywych.
rf_le_ale = rf_le_explainer.model_profile(N=1000, random_state=42, type='accumulated')
rf_ohe_ale = rf_ohe_explainer.model_profile(N=1000, random_state=42, type='accumulated')
xgb_le_ale = xgb_le_explainer.model_profile(N=1000, random_state=42, type='accumulated')
gb_le_ale = gb_le_explainer.model_profile(N=1000, random_state=42, type='accumulated')
rf_le_ale.plot([rf_ohe_ale, xgb_le_ale, gb_le_ale])
W przypadku wykresów ALE krzywe dla poszczególnych zmiennych i modeli są podobne do wcześniej otrzymanych profili PD, podobnie można interpretować wpływ wartości poszczególnych zmiennych objaśniających.
Krzywe ALE są również podobne między sobą, zatem możemy wnioskować o zgodności modeli, przede wszystkim obu modeli lasu losowego.
Wygenerujmy jeszcze porównanie wykresów ALE dla naszego głównego modelu z wcześniej otrzymanymi i przeanalizowanymi wykresami PDP.
rf_le_pdp.result._label_ = 'PDP'
rf_le_ale.result._label_ = 'ALE'
rf_le_ale.plot(rf_le_pdp)
W przypadku części zmiennych krzywe wygenerowane przy użyciu obu metod niemalże pokrywają się. Natomiast są też zmienne, gdzie widać różnice w wartościach predykcji (największe dla lead_time i country, a więc jednych z najważniejszych zmiennych), ale kształt krzywej pozostaje ten sam. Równoległość profili sugeruje i pozwala nam wnioskować, że wykorzystywany model jest addytywny ze względu na te zmienne wyjaśniające (na podstawie EMA).
Jednocześnie potwierdzają się wcześniej opisane wnioski dotyczące wpływu poszczególnych zmiennych.